Vaccination intervention on epidemic dynamics in networks 
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Vaccination is an important measure available for preventing or reducing the spread of infectious diseases. 
In this paper, an epidemic model including susceptible, infected and imperfectly vaccinated compartments is 
studied on Watts-Strogatz small-world, Barabasi- Albert scale-free, and random scale-free networks. The epi- 
demic threshold and prevalence are analyzed. For small-world networks, the effective vaccination intervention 
is suggested and its influence on the threshold and prevalence is analyzed. For scale-free networks, the thresh- 
old is found to be strongly dependent both on the effective vaccination rate and on the connectivity distribution. 
Moreover, so long as vaccination is effective, it can linearly decrease the epidemic prevalence in small-world 
networks, whereas for scale-free networks, it acts exponentially. These results can help adopting pragmatic 
treatment upon diseases in structured populations. 

PACS numbers: 89.75.Hc, 87.23.Ge, 87.19.X- 



I. INTRODUCTION 

Mathematical characterization of infectious diseases has 
contributed greatly to getting insight on transmission pat- 
terns of a disease in host populations, as well as on pub- 
lic health policies to prevent, reduce, and possibly erad- 
icate the disease [1-3]. Classical epidemic models usu- 
ally assume that either individuals do not have immu- 
nity to infection (the susceptible-infected-susceptible (SIS) 
model) or experiencing infection with permanent or tem- 
porary protection against it (susceptible-infected-recovered 
(SIR) and susceptible-infected-recovered-susceptible (SIRS) 
models). However, there is increasing evidence that most in- 
fections, such as pertussis and tuberculosis, can provide only 
partial immunity and spread among seropositive individuals, 
regardless of a reduced transmission rate. In view of this fact, 
vaccination was introduced into mathematical compartmental 
models which is often represented by a transfer between the 
susceptible and removed classes [4—9]. Whether vaccination 
is inoculation or education, typically it reaches only a fraction 
of the susceptible populations and is not perfectly effective. 
Thus, a backward transfer must be considered because vacci- 
nated individuals may return to be susceptible or become di- 
rectly infected. When these aspects are included in the model, 
rich dynamical behaviors may arise, such as backward bifur- 
cation and bistability [5, 6]. 

Previous studies of mathematical models incorporating vac- 
cination either ignore the population structure or treat popula- 
tions as distributed in a regular medium, that is, all the in- 
dividuals have the same probability of contacting the others. 
Recently, classical epidemic models have been extended in 
many ways (e.g., to study the disease spreading in a popula- 
tion divided into subgroups which may influence each other 
[10]). Especially, a great source of inspiration to mathemat- 
ical epidemiology has been provided by the network theory 
whose nodes represent individuals and links stand for inter- 
actions among them [11-14]. The structure of the underlying 
network (e.g., the degree distribution) may strongly influence 
spreading dynamics [15-27]. For instance, in scale-free (SF) 



networks, characterized by degree distributions with power- 
law behavior P(k) ~ fc~ 7 , the statistical relevance of hubs 
makes the network highly permeable to disease propagating 
[16]. This radical change in the behavior of the processes sug- 
gests that the standard epidemiological frameworks should be 
carefully revisited. 

The mathematical compartmental theory focuses on epi- 
demic equilibria and their stability. The network-based mod- 
eling, however, pays much attention to the underlying contact 
structure among individuals. The goal of this paper is to in- 
vestigate the influence of vaccination on disease spreading. 
Different from the classic study of the SIS model with vacci- 
nation by the compartmental theory which focuses on the sta- 
bility of equilibria [5], the present work revisits the model on 
Watts-Strogatz (WS), Barabasi-Albert (BA), and random SF 
networks concentrating on the epidemic threshold and preva- 
lence. The choice of this model is based on three factors as 
follows, (i) The SIS epidemic framework has been widely 
used in modeling disease spreading within a population. Each 
individual is simply assumed to have only one of the two 
states: susceptible (S) and infected (I). Each susceptible in- 
dividual gets infection with a transmission rate a once it con- 
tacts an infected one. Meanwhile, infected individuals recover 
and become susceptible again with a recovery rate (3. Then 
the process of disease transmission flows as S— ^1— ^S. (ii) Im- 
munization of population through vaccination strategy is an 
important and feasible practice with obvious implications for 
the public health. At the population level, it is interesting to 
determine the critical vaccination rate necessary for eradicat- 
ing diseases or preventing infection, and to investigate how 
vaccination affects the epidemic prevalence in the steady state 
on different networks. In this paper, to study possible effects 
of vaccination on epidemic dynamics in different networked 
populations, a vaccinated (V) state is introduced into the SIS 
model by vaccinating the susceptible individuals with a vac- 
cination rate tp, corresponding to the transition S— >V. (iii) In 
the real world, there are various types of vaccines: some may 
offer temporary immunity; vaccines may not possess 100% 
efficacy (leaky vaccines) [28], and finally, on most occasions, 
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FIG. 1: Flowchart of the SIS model with vaccination. Susceptible 
nodes are infected by their infected neighbors at a per capita rate 
a and are vaccinated at a per capita rate <p. Infected nodes recover 
to be susceptible at a per capita rate j3. Vaccinated nodes become 
susceptible at a per capita rate <j> and are infected at a per capita rate 
8a due to imperfect vaccination. 

vaccination may cover only a fraction of susceptible individu- 
als. Therefore, the vaccinated individuals return to the suscep- 
tible class with a resusceptibility rate <f> as the vaccine wears 
off, or directly get infected with a reduced transmission rate 
5a, where 5 denotes the degree to which the vaccine-induced 
protection against infection is inefficient. Thus, the vaccinated 
class flows in two directions S<— V— >I. 

To study the SIS model with vaccination on networks, both 
analytical calculations and numerical simulations are carried 
out. Depending on network structures, two types of epi- 
demic thresholds and corresponding prevalence behavior are 
obtained. For WS networks, there is a nonzero threshold 
similar to that obtained in the classic compartmental model. 
Whereas for SF networks, the threshold is quite different 
which is strongly related both to the vaccination rate and to 
the node-degree distribution. While the disease is endemic 
in the network, it is found that vaccination intervention con- 
tributes to linearly reducing the prevalence in WS networks, 
whereas in SF networks it functions exponentially. 

The rest of the paper is arranged as follows. In Sec. II the 
SIS model with vaccination is introduced. Then the model is 
studied on WS, BA, and random SF networks in Sees. Ill, IV, 
and V, respectively. Finally, conclusions are given in Sec. VI. 

II. THE SIS MODEL WITH VACCINATION 

Vaccination scheme has always been a very important and 
effective way for preventing or controlling infectious diseases. 
In reality, vaccines hardly cover the whole population, and 
only remain effective for a finite period of time, and are diffi- 
cult to guarantee a perfect protection from infection. Taking 
all these points into consideration, Kribs-Zaleta and Velasco- 
Hernandez introduced a vaccinated state into the SIS model 
to theoretically study the possible effects of vaccination on 
epidemics [5]. The mathematical compartmental model in 
Ref. [5] neglects the population structure and assumes all indi- 



viduals have the same contact rate. Therefore, it is intriguing 
to inspect the vaccination program for different network struc- 
tures. 

In this paper the SIS model with vaccination is formulated 
on the static network framework, where nodes represent in- 
dividuals and links stand for the contacts among individuals 
along which a disease can spread. At each time step, each 
node exists in only one of the three states: susceptible, in- 
fected and vaccinated. A disease spreads in the network fol- 
lowing mechanisms below (as shown in Fig. 1): a susceptible 
node will be infected with the transmission rate a once it is 
connected to an infected one; an infected node is cured and 
become susceptible again with the recovery rate j3; according 
to the random vaccination strategy, each susceptible node gets 
vaccinated with the vaccination rate ip and each vaccinated 
node returns to the susceptible class at the resusceptibility rate 
<j) as the vaccine wears off; due to the imperfect immunity, 
a vaccinated node will be infected with a reduced infection 
rate 5a, where the parameter 5 measures the inefficacy of the 
vaccine-induced protection against infection ((5 = and 1 re- 
spectively represent completely effective and utterly invalid, 
but < 5 -c 1 holds for most cases [28]). In present work it 
is always assumed that 5 is sufficiently small. The epidemic 
dynamics is determined by five parameters, a, (3, ip, <j>, and 5. 
For convenience, two ratios A, ?y are particularly denoted as 
A = a//3 and r\ — (p/<f) [29]. 

In the theoretical study of epidemiology, there are two 
paramount indicators. One is to formulate the epidemic 
threshold, which determines whether the infection breaks out 
in the population and results in an EE or dies out eventually 
corresponding to a DFE [30]. The other is to predict the epi- 
demic prevalence. The present work make a comprehensive 
study of the SIS model with vaccination on WS, BA, and ran- 
dom SF networks, employing the mean-field (MF) approach 
and computational simulations. As will be seen below, vacci- 
nation indeed has a great influence on the epidemic dynamics 
over such networks. 



III. THE MODEL ON WS NETWORKS 

The WS network [31], as a reference of homogeneous net- 
works, can be constructed as follows. Start from a ring of 
N nodes, where each node is connected symmetrically with 
its 2K nearest neighbors. Then, every link connected to a 
clockwise neighbor is rewired to a randomly chosen node with 
probability p. After the whole sweep, a WS network with the 
average connectivity (k) = 2K is generated. 

The WS network is a typical example of networks charac- 
terized by a narrow degree distribution, in which each node's 
degree closes to (k). Let s(t), pit), and v(t) be the densities 
of susceptible, infected and vaccinated individuals at time t, 
respectively. Obviously, they satisfy the normalization con- 
dition s(t) + p{t) + v(t) = 1. Therefore, a set of coupled 
differential equations can be established following the MF ap- 
proach [17]: 
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-p(t) = -/3p(t) + a{k)s(t)p{t) + Sa(k) [1 - p(t) - s(t)}p(t), 
^s(t) = pp{t) - a{k)s(t)p(t) + <(>[1 - p(t) - s(t)] - ips{t). 

I 



(la) 
(lb) 



The first term on the right-hand side (rhs) in Eq. (la) accounts 
for the recovery process from the infected class, which is pro- 
portional to the recovery rate j3, the average density p(t) of in- 
fected nodes. The second term on the rhs in Eq. (la) denotes 
the newly infected nodes transferred from susceptible ones. 
It is proportional to the density s(t) of susceptible nodes, the 
transmission rate a, the average number of neighbors (k), and 
the probability p(t) that a randomly chosen neighbor is in- 
fected. Similarly, the third term on the rhs in Eq. (la) consid- 
ers the probability that a node is vaccinated [1 — p(t) — s(t)], 
and gets infection. The probability of this last process is pro- 
portional to the vaccine-reduced transmission rate 8a, the av- 
erage number of neighbors (fc), and the probability p{t) that a 
randomly chosen neighbor is infected. On the rhs in Eq. (lb), 
the third term accounts for the increment in the susceptible 
class result from the transition V^S, which is proportional to 
the resusceptibility rate 0; and the fourth term accounts for 
the probability of the vaccination process S— >V, which is pro- 
portional to the vaccination rate <p. 

It should be stressed that the MF approach (Eqs. (la) and 
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(lb)) is equivalent to the mass-action law (system (1) in 
Ref. [5]) with the adaptation of the reaction rates to include 
the average connectivity (k). Thus, one can obtain similar re- 
sults for the epidemic threshold and equilibrium stability. It 
is due to the homogeneity of the WS network, which is also 
the case for the SIS model in Ref. [17]. In accordance with 
the results in Ref. [5], there is either a forward bifurcation or 
a backward one, depending on the epidemiological parame- 
ters: (1) in case of the model exhibiting the forward bifurca- 
tion, there is only one globally stable EE as the transmission 
rate is above the epidemic threshold, below which the DFE 
is the only attractor; (2) in case of the presence of the back- 
ward bifurcation, there are multiple endemic equilibria (MEE) 
- meaning that there are two or more EEs in the steady state - 
existing between a sub-threshold and the epidemic threshold, 
meanwhile both the DFE and the lower EE are locally stable. 

In the following much attention will be paid on the epi- 
demic prevalence in the steady state. Imposing the stationary 
conditions j^s(t) = and j^p(t) = yields 



(1- (3)p + 8a(k)(l- p)p + a(k)(l-6) 



Qg - 4>)p 2 + 4>p 

a(k)p + ip + <f) 



f(p) 



(2) 



r 



for density of infected nodes in steady states. Notice that 
/(0) = and f(l) < 1, Eq. (2) has a nonzero solution on 



p=0 



> 1, which defines the 
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the interval (0, 1) only if f'(p) 
epidemic threshold 

_ ^ + 1 

c &ip + (j> (k) Sr] + 1 (k) ' 
On the other hand, Eq. (2) can be rewritten as 

F(p) = Ap 2 + Bp + C = 0, 
with coefficients 

A = 5a(k), 

B = 5^ + <j) + Sp-5a(k), 

c = -[(Sip + <!>)- + 



The solutions to Eq. (4) correspond to equilibria of system (1) 
for given A. 



(i) A > A c . In this case, C < always holds. Since F(0) = 
C < and F(l) = A + B + C = 6f3 + (ip + (f>)/(X{k)) > 0, 
system (1) has a unique EE, 



P 



v*B 2 - AAC - B 
2A 



(5) 



In contrast to special solutions obtained in Ref. [5], this ex- 
pression is general. As A — > A c , C closes to 0. From Eq. (5), 
it follows that p — > 0. Ignoring the second order term of p in 
Eq. (4), one has 



P 



dtp + cj) 



X-K 



5ip + 4> + 5(3 - 5a(k) A c 



(A -A c ). (6) 



In fact, as given in Appendix A, applying Taylor series ex- 
pansion to the square root part of the first term at 6 = and 
omitting the higher order correction in S, Eq. (5) can be sim- 
plified to 
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FIG. 2: (Color online) Bifurcation diagrams of p as a function of A in 
the WS network with the average connectivity (k) = 10 for various 
S: 0.001 (a), 0.01 (b), and 0.1 (c). In each chart, black and red curves 
respectively represent stable and unstable branches of system (1). 
Parameter values: /3 = 0.002, ip = 0.001, and <j> = 0.0002. 
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FIG. 3: (Color online) Illustration of the dynamic behavior on the 
8-r) plane with parameters f3 = 0.002 and <j> — 0.0002, showing 
the necessary condition for bistability in the WS network with (k) = 
10. The white (grey) region corresponds to the single stable state 
(bistable states). The blue (red) line is upper (lower) bound of rj for 
given S. Only if r\ is inside the range between the upper and lower 
bounds can the MEE occur. Otherwise, there is only one attractor 
(either the DFE or the EE). 



(ii) A < A c . In this case, it is impossible to obtain the 
general solution of the epidemic prevalence. Following the 
parametric analysis in Ref. [5], one finds that there are two 
different EEs in the regime Ab < A < A c on the premise of 
A a < Ab, where A a = S ^^^g^ corresponds to the condition 



B = 0, and A b 



l 

<fe) 



responds to B 2 
the condition A a 



1 



5/3 



- AAC = 0. This finding suggests that under 
< Ab, which is equivalent to 



(Sip + <j>) 2 < (35(1 - 6)<p, 



(8) 



there emerges a sub-threshold Ab (the persistence threshold, 
above which an already established epidemic can persist [32]) 
and an epidemic threshold A c (the invasion threshold, which 
still denotes the critical parameter value for invasion of new 
diseases). In the bifurcation diagram, this sub-threshold corre- 
sponds to a saddle-node bifurcation, and the epidemic thresh- 
old corresponds to a forward bifurcation, and consequently 
these two thresholds together exhibit first order transitions be- 
tween the healthy phase (without disease) and the endemic 
phase (with disease) (see Fig. 2(c)). This reveals the hys- 
teresis effect caused by the introduction of vaccination into 
the infectious disease. However, in other cases, there is only 
the invasion threshold, and hence only the forward bifurcation 
(see Figs. 2(a) and 2(b)). 

Clearly, whether the interval (Ab, A c ) is a bistable region or 
not in the bifurcation diagram of p as a function of A is com- 
pletely determined by the condition (8) which can be rewritten 
as 



uV + [2-^(1 -S)]Sri- 



1 < 0. 



(9) 



Only if 4(1 — 8) > 4 can the inequality has solutions on 
the interval (771,772) C (0,oo), which indicates that S max = 
1 — 4^. On the contrary, there is only one single stable 
state if /3/(f> < 4. Figure 3 gives an illustration. In case of 
(3 = 0.002 and = 0.0002, the ratio is /?/</> = 10. To en- 
sure 4(1 — 5) > 4, it demands S < 0.6. In particular, at 
S = 0.001, the upper and lower bounds of 77 are 127.18(1) and 
7892.81(6), respectively. On the other side, given 77 = 5.0, 
only if 0.02(7) < 6 < 0.51(1) can system (1) experience the 
bistable states. From calculation in Appendix B, r\ reaches the 
minimum 2/3 at <5 = 3/8. The presence of such a bistable 
region highlights an important but unexpected influence of 
vaccination on disease spreading. The bifurcation diagrams 
in Fig. 2 correspond to three different values of 6 (= 0.001, 
0.01, and 0.1) for 77 = 5.0. At S = 0.1 the system exhibits 
the bistable phenomenon, and in order to wipe out the disease, 
one must ensure that A < 0.37(4) rather than A < 0.4. 

The emergence of MEE gives rise to complexity in the vac- 
cination intervention. In the real world, it is interesting to 
study the effective vaccination and its influence on the epi- 
demic prevalence. To do that, one can consider system (1) 
with proper choice of (3, <fi, 77, and 6, respectively, ensuring 
that only the forward bifurcation occurs. In the following, S 
is fixed at a relatively small value S = 0.001 with (3/<j> = 10 
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TABLE I: Simulation values of Do, Di calculated by applying the Levenberg-Marquardt algorithm [33] to the least squares curve fitting on 
the simulation data plotted in Figs. 4(b) and 4(c), with the general function in the form of p = Do — Di 2 ^-. The quantitative comparison is 
also demonstrated, which shows a good agreement between the numerical simulation and the analytical prediction by Eq. (7). 
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FIG. 4: (Color online) Infected densities p in the WS networks as 
functions of A (a), 1/A (b), and r\ (c), respectively. Solid lines are 
analytical solutions to Eq. (5). Dash lines are theoretical prediction 
by Eq. (7). Parameters values: /3 = 0.002, <f> = 0.0002, and S = 
0.001. 



and r\ < 127.18(1), where only a globally stable EE exists if 
A > A c or a globally stable DFE arises if A < A c . Accord- 
ing to the model definition, A and 77 are comparably important 
parameters which affect the global spread of the infection. In 
Fig. 4 both analytical and numerical results of p as a func- 
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FIG. 5: (Color online) 1/A C vs 1/(77 + 1) in the WS networks for dif- 
ferent connectivities. Solid lines correspond to solutions of Eq. (3). 
Parameter values: /3 = 0.002, <j> = 0.0002, and S = 0.001. 



tion of A and rj in the WS network are present, respectively. 
Simulation of the SIS model with vaccination on the WS net- 
work is carried out with parameters = 10 5 and p = 0.1. 
The fraction of initial infectious seeds is 0.1% and the preva- 
lence p in the steady state is averaged over 10 different real- 
izations of the model on each of 10 different initial network 
configurations. Each realization goes through 2 x 10 4 time 
steps. The thresholds A c in Fig. 4(a) are 0.77(4), 0.62(5), 
and 0.52(3), corresponding to the average degrees (k) = 8, 
10, and 12, respectively, which agrees with the prediction of 
Eq. (3). Moreover, linear behaviors are shown from both the 
simulation results and the theoretical predictions in Figs. 4(b) 
and 4(c). 

To examine the accuracy of analytical prediction by Eq. (7), 
a quantitative comparison is made in Table I, where the 
Levenberg-Marquardt algorithm is applied to the least squares 
curve fitting on the simulation data, with the general function 
in the form of p = Dq — D\ (77 + 1)/ A, where Do and D\ are 
positive constants. In Table I, the numerical value Dq ranges 
from 1.00(3) to 1.07(2), matching the theoretical prediction 
Dq = 1; and the numerical Di is also in good agreement with 
the prediction D\ = Thus, given all the epidemiolog- 
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ical parameters, the larger the average node-degree (k), the 
harder for disease to break out, and if it prevails, the higher the 
level of infection p forms (as shown in Fig. 4). On the other 
hand, in a fixed WS network, as the ratio (77 + 1)/A increases, 
p linearly diminishes. This implies that the competition be- 
tween the transmission process and the vaccination campaign 
leads to a linear decrease of the prevalence in networks with a 
narrow degree distribution, as shown in Fig. 4(c). Notice that 
for sufficiently small 5, as assumed in present work, Eq. (7) 
turns into p as 1 — A c /A. In form, this scaling behavior of 
p resembles that in the SIS model without vaccination [17], 
where the threshold is A c = l/(k) . This similarity, how- 
ever, reveals that the linear effect of vaccination on the preva- 
lence is in essence due to the fact that vaccination program 
increases the epidemic threshold by 77 times. To get further 
information, the inverse of A c as a function of the inverse of 
r/ + 1 is depicted in Fig. 5. Since Eq. (3) can be rewritten 
as 1/A C = — fy/iv + 1) + S\> simulation results ver- 

ify this linearity. Hence, the more effectively the vaccination 
intervenes on the disease, the more difficultly it outbreaks. 



IV. THE MODEL ON BA NETWORKS 



The BA network [34], as a prototype of heterogeneous net- 
works, can be built as follows. Start from a a set of mo nodes, 
which are completely connected. At each time step, a new 
node is added to the existing network, bringing m(< mo) 
new links connecting to old nodes with degree preference. Af- 
ter iterating this procedure a sufficient number of times, a B A 
network is obtained, consisting of TV nodes with the node- 
degree distribution P(k) = 2m 2 k~ 3 and the mean node- 
degree (k) = 2m. 

The heterogeneity of the connectivity distribution inherent 
to BA networks induces strong fluctuations, so systems (1) 
should be modified accordingly. Denoting by s k (t), Pk{t), 
and Vk{t) the relative densities of susceptible, infected and 
vaccinated nodes with degree k at time t, respectively, which 
satisfy the normalization condition s k (t) + Pk (t) + Vk (t) = 1, 
the MF equations now read as: 



g- t Pk(t) = ~Ppk{t) + aks k (t)Q{p(t)) + 5ak[l - p k (t) - s k (t)]Q(p(t)), 
d 

g- t s k{t) = Ppkif) - aks k (t)Q(p(t)) + 4>[1- p k {t) - s k {t)] - ips k (t). 



(10a) 
(10b) 



The first term on the rhs in Eq. (10a) considers that a node 
of degree k is in the infected state with probability p k (t) and 
recovers from infection at the recovery rate f3. The second 
term on the rhs in Eq. (10a) considers the probability that a 
node with k links is in the susceptible state s k (t) and gets in- 
fection via a neighbor. The probability of this last event is 
proportional to the transmission rate a, the number of neigh- 
bors k, and the probability 8 (/?(£)) that any given link points 
to an infected node. Similarly, the third term on the rhs in 
Eq. (10a) considers that a node with k neighbors is in the vac- 
cinated state with probability [1 — p k (t) — s k (t)] and gets 
infection via a neighbor at the vaccine-reduced transmission 
rate da. On the rhs in Eq. (10b), the third term considers that 
a node with degree k is in the vaccinated state with probability 
[1 — p k (t) — s k (t)] and returns to the susceptible class at the re- 
susceptibility rate <f>; and the fourth term considers that a node 
of degree k is susceptible with probability s k (t) and gets vac- 
cinated at the vaccination rate (p. For uncorrected networks, 
the probability 9 is [17] 



kp(k) 



Pk- 



(11) 



Since SF networks have no correlations under the constraint 
that the maximum possible degree has a cutoff scaling at most 
as k c (N) ~ TV 1 / 2 [35]. In order to ensure an uncorrelated BA 
network, this restriction on the maximum degree is imposed in 
present work. Imposing the stationary conditions ^p k {t) = 



and 4is k (t) = yields 



Pk 



ak<d(5ip + <f> + 5akQ) 



vkQ(5(p + 4> + SakQ) + f}(y> + cf> + SakQ) 



(12) 



Combining Eqs. (11) and (12), one obtains a self-consistency 
equation, 



e 



1 

W) 



5(e). 



kP(k)ak@(5(p + (f> + Sake) 
vkQ(5ip + + SakQ) + p(ip + <f) + SakQ) 

(13) 



Obviously, there is a trivial solution O = which leads to 
p = 0. Notice that not only g(0) = and g(l) < 1, but 
also g'(Q) > 0, and ,g"(6) < in the limit 5^-0, only if 

<?'(0) > 1 can Eq. (13) have a nontrivial solution on the 

interval (0, 1), which yields 



Xc = 



tp + <f> (k) _ T} + 1 (k) 



5if + (j){k 2 ) 5r]+l(k 2 )' 



(14) 



In infinite-sized BA networks, the second moment of the con- 
nectivity distribution is unbounded, i.e., (k 2 ) — > 00, which 
induces A c = 0. So the infection can always prevail among 
the population, no matter what the effective transmission rate 
is. Whereas for finite-sized BA networks, there exists a maxi- 
mum degree k c , which controls the bound of the connectivity 
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fluctuations, inducing a nonzero threshold [15, 18]. From now 
on, the size of the BA networks is assumed to be finite, and all 
the possible values of node degrees are k = m, m + 1, . . . , k c . 
By computing the Jacobian matrix of the DFE {(p^ = 



0, Sk 



<p+<j> ■ 



)}k=m °f svstem (10), one finds that the basic 



reproductive number [1] is Rq = A Kj-J- = A/A c , which 
denotes the expected number of secondary infections caused 
by a single infected individual in a completely susceptible 
population. Accordingly, the DFE is locally asymptotically 
stable if A < A c , while unstable if A > A c meaning invasion 
is always possible. As long as A > A c there exists a positive 
solution p E (0, 1) corresponding to the EE which is locally 
asymptotically stable. 

Since 9 approaches as A closes to A c , and by neglecting 
all higher order corrections in 8, Eq. (12) is in form analo- 
gous to Eq. (8) in Ref. [17], and hence one expects the similar 
critical behavior given by Eq. (14). Compared with the SIS 
model on BA networks [16, 17], the presence of vaccination 
has the effect of multiplying the epidemic threshold by a factor 
(77 + 1)/ (Srj + 1), i.e., enlarging by nearly 77 times (as S — > 0). 
This suggests that vaccination might play a significant role in 
preventing or reducing the infectious disease. The greater the 
vaccination rate is, the bigger the epidemic threshold is, and 
hence the harder the disease erupts. 

Neglecting the second order term in 9, Eq. (12) can be sim- 
plified as 



Pk 



\Qk + T 1 + l' 



(15) 



Given the epidemiological parameters, a nodes with higher 
degree is more likely to get infected. Substituting this ex- 
pression into Eq. (13) and treating k as a continuous variable 
yields 



9 « toA9 



1 



k \Ok 



dk 



V 



(16) 



which gives rise to the solution 



QRj(??+ l)e-(^ _ {n+l)l ^ 
Am 

Finally, at lowest order in A, the epidemic prevalence related 
to the EE is 



)/\r, 



(18) 



Computational simulations for the epidemic model are per- 
formed on the BA networks with the network size N = 10 s . 
Each of the simulation data is obtained by averaging over 10 
different realizations of the model on each of 10 different net- 
work configurations. Each realization goes through 2 x 10 4 
time steps. As shown in Fig. 6, there is a deviation between the 
simulation results and the analytical calculations, especially in 
the large prevalence regime. It is due to the fact that Eq. (12) 
is simplified by neglecting the highest order in 9. As p is 
relatively large, 9 2 is actually not negligible. Despite this, 




FIG. 6: (Color online) Densities of infectious nodes p in the BA 
networks: (a) as a function of 1/A for 77 = 5.0 and various (fc); 
(b) as a function of 77 for A = 0.4 and various m. Solid lines are 
theoretical prediction by Eq. (18). The other parameters are j3 = 
0.002, (j> = 0.0002, and 5 = 0.001. 
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TABLE II: Simulation value of Eo calculated by applying the 
Levenberg-Marquardt algorithm to the least squares curve fitting on 
the simulation data plotted in Figs. 6(a) and 6(b). The quantitative 
comparison is also demonstrated, between the fitting value Eo and 
the analytical prediction 1/m by Eq. (7). 



the simulation support the calculation by the same exponen- 
tial decaying in the scaling behavior, i.e., p ~ e - E o(v+i)/^^ 
where Eo is a constant. The numerical comparison is also 
made between the fitting value Eq and the analytical estima- 
tion 1 jm in Table II, showing a relatively small variance. For 
BA networks, both in simulation and in theory, the prevalence 
decays exponentially, i.e., p ~ e ^i r i+ 1 )/ Xm ^ Vaccination has 
an effect of accelerating by 77 times the exponential decreasing 
of the prevalence. This finding suggests that the vaccination 
intervention on a disease can efficiently reduce an endemic to 
a lower level, though the heterogeneity in degree distribution 
causes a vulnerability to disease outbreak in BA networks. 
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V. THE MODEL ON RANDOM SF NETWORKS 

In this section, the analysis for the SIS model with vac- 
cination on BA networks will be generalized to random SF 
networks with arbitrary exponent 7 > 2. Following the idea 
proposed by Newman et al. [36], the random SF networks 
can be generated as below. First, a priori random integers se- 
quence, each of which represents the degree of a node, drawn 
from a normalized distribution 



if k < k c . 
otherwise, 



(19) 



where m and k c are respectively assumed to be the minimum 
and the maximum values of the degree among all the nodes, 
and k c 3> m. Notice that in order to get uncorrected ran- 
dom SF networks, the restriction on the maximum degree [35] 
k c {N) ~ N 1 / 2 is imposed. Then, node i with degree fej is 
picked out randomly from the sequence and connected to oth- 
ers until its degree quota hi is realized. Duplicate connections 
are avoided. This process is repeated throughout all the ele- 
ments of the sequence, and finally a network is chosen uni- 
formly at random from the set of all graphs with that degree 
sequence. Assuming k changes continuously and the average 
connectivity is thus 



kP(k)dk « - l -m. 



(20) 



For any connectivity distribution in random SF networks, 
one can employ directly the analytical treatment in BA net- 
works. That is to say, the MF results in the BA network are 
applicable to the random SF network. According to Eq. (14), 
the epidemic threshold is zero if 7 < 3 in the thermodynamic 
limit for random SF networks. Whereas for 7 > 3, substitut- 
ing Eq. (20) into Eq. (14) yields a non-zero threshold 



Ac 



(77 + 1X7-3) 
m(5ri + 1)(7 - 2) ' 



(21) 



(l-3)(7-3) 



< for any 6 G (0, 1), en- 



Since J^A C — ( 5 , )+1 )2 (7 _ 2 ) rn 
hancing the effective vaccination rate can prevent epidemics 
from spreading through the population. For any exponent 7 in 
random SF networks, combining Eqs. (13) and (20), one has 



e 



A'6(7-l)m^- 1 f fcc fc 2 -T 



(k) 



A'Gfc + 1 



with 



A' 



V 



1 



dk, (22) 



(23) 



Due to the existence of the parameter 7, it is difficult to obtain 
the explicit solution of Eq. (22). However, one can roughly 
estimates and p using the first mean value theorem. In this 
way, one has 



6 = 



A'6(7-2)m 
1 + A'Gfii 



7-2 



k^dk, 



(24) 



where fli is a finite constant, m < Qi < k c . Thus, the solu- 
tion is 



e 



1 



A'(7-2)m 7 ~ 2 fcg~"' 



1 

A'Qi 



3-7 
A'(7-2)m 
7-3 



- 1 



if 2 < 7 < 3, 



if 



7 > 3. 



(25) 



The prevalence p can also be written as 

p= [ kC p(k)p k dk* "^'^ (l- 
Jm r2 2 ( 7 -2) 



1 



1 + 



), (26) 



where 17 2 is a finite constant, m < 0,2 < k c . According to 
Eq. (25), the behavior of p depends on 7. 

(i) 2 < 7 < 3. In this case, for any A > 0, A'0O 2 > 1. 
Implementing logarithm operation on Eq. (26) yields 



In p w In 



771(7 ~~ 1) 



1 



n 2 (7-2) l + A'6fV 
Combining this with Eq. (23), one has 

where V\ is a constant, defined by 

fii(3-7) 



f2 2 mT- 2 ^- 7 ( 7 -2)' 



(27) 



(28) 



(29) 



(ii) 7 > 3 while 7 ^> 3. According to Eqs. (21) and (25), 
for any A A c , X'QQ.^ 0. One can obtain the prevalence 
similar to case (i) 



P 



-U 2 (77+1)/A 



where the coefficient v-i reads as 



''2 



Oi(7-3) 
n 3 m(7-2)' 



(30) 



(31) 



(iii) 7 3> 3. In this case, the connectivity distribution de- 
cays so fast that it tends to a homogeneous networks. One 
would expect to obtain the similar qualitative behavior as in 
Sec. III. 

Simulations of the SIS model with vaccination on random 
SF networks are performed to compare with the theoretical 
analysis. The simulated networks range from N = 10 5 to 
N = 10 6 and the minimal degree of nodes is m = 5. Figure 7 
shows the epidemic threshold A c as a function of the alge- 
braic expression (7 — 3)/(7— 2). Closed squares represent 
numerical data and the solid line corresponds to the predic- 
tion of Eq. (21). One notices the good agreement between the 
computer simulation and the analytical calculation. Figure 8 
depicts the behaviors of p as a function of 1 /A (Fig. 8(a)) and 
as a function of 77 (Fig. 8(b)), respectively. It is clear that ei- 
ther for A > (the case of 2 < 7 < 3) or for A > A c (the 
case of 7 > 3 and 7 ^> 3), the stationary density p of infected 
nodes in the random SF networks decays exponentially, i.e., 
p ~ e ~ v ( r i+ 1 )/ x ^ where v is a positive constant, which is de- 
termined by Eq. (29) or Eq. (31). On the contrary, at 7 — 40, 
as shown in the inset of Fig. 8, p decreases linearly similar to 
the behavior observed from WS networks. 
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FIG. 7: (Color online) Effective transmission threshold A c as a func- 
tion of (7 — 3)/ (7 — 2) in the random SF network. The full line 
corresponds to the analytical calculation of Eq. (21). Parameter val- 
ues: 77 = 5.0, f3 = 0.002, (j) = 0.0002, and r5 = 0.001. 
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FIG. 8: (Color online) Semi-log plots of the persistence p in random 
SF networks as a function of 1/A (with 77 = 5.0) (a) and 77 (with 
A = 0.5) (b) for various values of 7: 2.5, 3.5, and 4.5 (from top 
to bottom). The insets of (a) and (b) respectively display the linear 
dependence of p as a function of 1/A and 77 for 7 = 40. Parameter 
values: /3 = 0.002, = 0.0002, and 6 = 0.001. 



a networked SIS model with vaccination, where vaccines that 
attempt to reduce susceptibility to infection is characterized 
by three parameters in the model: coverage (represented by 
tp), waning period (represented by (j>), and efficacy (repre- 
sented by 5). Since S is intrinsically related to the quality 
of the vaccine, much attention has been paid to the parame- 
ter 77 (the ratio of ip to <f>) in the vaccination intervention on 
infectious diseases, as well as the role of the ratio A (the ra- 
tion of a to f3) in epidemic spreading. With the frameworks of 
the MF approach and elementary means, the model has been 
studied on WS, BA, and random SF networks. The analysis 
of thresholds and prevalence demonstrated the significant ef- 
fects of the vaccination on the epidemic dynamics as well as 
the structures of the underlying networks. 

In the WS networks, since the MF model is equivalent to 
the classic compartmental model in Ref. [5] with the adapta- 
tion of reaction rates by the average connectivity, the threshold 
behavior and equilibrium stability are similar to the literature. 
The threshold A c is defined by Eq. (3), above which there is 
only one globally stable EE and below which the model may 
exhibit MEE for certain epidemiological parameters. As to 
the prevalence, rather than special solutions obtained in the 
compartmental model, this paper gives the general one for 
the steady endemic state, which scales as p ~ — (77 + 1)/A. 
Thus, the effective vaccination can linearly decrease the en- 
demic level in homogeneous networks, although vaccination 
intervention may give rise to the backward bifurcation in these 
networks. 

In the SF networks, however, the system shows very differ- 
ent behavior. The threshold A c is defined by Eq. (14). Only 
for the SF network with the power-law distribution exponent 
2 < 7 < 3 in the thermodynamic limit can A c be zero. Oth- 
erwise, the system has a non-zero threshold for the SF net- 
works with any 7 > 3. In comparison with the WS network 
at the same average connectivity (k), A c in the SF network is 
smaller than that in the WS network. For any A > A c , the 
prevalence in the SF network scales as p ~ e - v i r i+ 1 )/ x , Thus, 
the vaccination can exponentially decrease the endemic level 
in heterogeneous networks. 

All these results are on the presumption that the underly- 
ing networks are static. For some diseases which spread too 
fast in comparison with change of the population structure, 
the present work may provide a preliminary theory for vac- 
cine control of infection. For other diseases, however, indi- 
vidual responses to infection plays an important role in either 
reducing the transmission rate or changing the contact struc- 
ture [37-39]. Hence, it is interesting to study vaccination in 
adaptive networks, which is left for future research. 



VI. CONCLUSION 

The study of vaccination in populations has to take into con- 
sideration not only vaccine-related parameters, but also social 
risk behaviors that may alter the expected predictions. To our 
knowledge, however, very few work addressed this problem. 
The present research integrated the both factors and studied 
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Appendix A: Simplification of Eq. (5) to Eq. (7) 

Since the inefficacy rate of vaccine 6 is assumed to be suf- 
ficiently small, it is possible to simplify the complex square 
root part in Eq. (5) via Taylor series expansion at the point 
5 = 0. First of all, rewrite Eq. (5) as 



where 



and 



Pi 



accompanied with 



P = Pi + P2, 

(k)a5 - (cp5 + 06 ■ 
2(k)a5 

92 2(k)aS 
h(6) = [h 1 (S)] 1/2 , 



(Al) 



(A2) 



(A3) 



(A4) 



hi(6) = (6y?+6+68-6a(k)) +48a(k)(6(p+d)-45/3(ip+6). 

(A5) 

So, rearranging each term, one has 



fn(5) = aS 2 + bS + 6 2 , 



(A6) 



where 



a= (tp + 8 - (k)a) + 4{k)aip, 
b = 2(k)a6 + 2ip6 - 2f36 - 4f3ip. 



Since the first order derivative of h(S) can be calculated as 

ti(S) = ^[h 1 (5)]- 1/2 (2aS + b), (A8) 



one has 



(k)a<b + ifid> - 86 - 20ip 
h'(0) = w HV tZ. (A9) 



Therefore, by employing Taylor series expansion at 5 = for 
h(8) with regard to 6, one gets 

h(6) = h(0) + h'(0)S + o(8 2 ) 

(k)a6 + ip6- 84>-2/3(p 



Combining Eq. (A10) with Eqs. (Al, A2, A3) gives rise to 

2(k)a6- 2/3(5- 2(3^5 



P ■ 



2(k)a6 



(All) 



which implies the simple relationship 

« 1 r, + l 
p « 1 



(A12) 



(k) A ■ 

Appendix B: Calculation of the minimum lower bound of 77 

Let x = 8/6 and substitute it into inequality (9), one ob- 
tains 



<5V + [2 - x(l - 6)}6r) + 1 < 0, 



(Bl) 



which has positive solutions if and only if 

{[x(l-5) -2] 2 -4}(5 2 > 0, <^> x(l - 6) > 4. (B2) 
The solutions of inequality (Bl) read 

m(S) <r)<ri 2 (6), (B3) 

where 



( x\ [x{1 - 6) - 2] - ^ [x{l 5) 2]2 - 4 ™ 



fXS W ~ 8) ~ 2] + ygj - S) - f - 4 

= Ys ■ (B5) 



The derivative of the lower bound rji (5) is 

d?7i 1 r x(x - 4) - x(x - 2)6 



dr/i 1 ( x[x — 



A /[ I (l-i)-2P-4 



(x-2)}. (B6) 



Let ^LL = 0, it follows that 



x[x-4-(x- 2)6} = {x- 2)vRl -5) - 2] 2 -4, (B7) 
which gives the extreme point 

d * = 27^) = 2if^y (B8) 

<5. (A10) Substituting 6* into Eq. (B4) yields the minimal lower bound 

4 4 

'Jlmin = Z 7 = fl 7- (B9) 
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